Jamming at Zero Temperature and Zero Applied Stress: the Epitome of Disorder 
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We have studied how 2- and 3- dimensional systems made up of particles interacting with finite 
range, repulsive potentials jam (i.e., develop a yield stress in a disordered state) at zero temperature 
and zero applied stress. At low packing fractions <f>, the system is not jammed and each particle can 
move without impediment from its neighbors. For each configuration, there is a unique jamming 
threshold, cj> c , at which particles can no longer avoid each other and the bulk and shear moduli 
simultaneously become non-zero. The distribution of <f) c values becomes narrower as the system 
size increases, so that essentially all configurations jam at the same packing fraction in the ther- 
modynamic limit. This packing fraction corresponds to the previously measured value for random 
close-packing. In fact, our results provide a well-defined meaning for "random close-packing" in 
terms of the fraction of all phase space with inherent structures that jam. The jamming thresh- 
old, Point J, occurring at zero temperature and applied stress and at the random close-packing 
density, has properties reminiscent of an ordinary critical point. As Point J is approached from 
higher packing fractions, power-law scaling is found for the divergence of the first peak in the pair 
correlation function and in the vanishing of the pressure, shear modulus, and excess number of 
overlapping neighbors. Moreover, near Point J, certain quantities no longer self-average, suggesting 
the existence of a length scale that diverges at J. However, Point J also differs from an ordinary 
critical point: the scaling exponents do not depend on dimension but do depend on the inter-particle 
potential. Finally, as Point J is approached from high packing fractions, the density of vibrational 
states develops a large excess of low-frequency modes. Indeed, at Point J, the density of states is a 
constant all the way down to zero frequency. All of these results suggest that Point J may control 
behavior in its vicinity - perhaps even at the glass transition. 

PACS numbers: 81.05.Rm, 82.70.-y, 83.80.Fg 



I. INTRODUCTION 



The nature of the glass transition has been called prob- 
ably "the deepest and most interesting unsolved problem 
in solid state theory." [| The nature of granular mate- 
rials has also been said to lead to equally deep ques- 
tions in statistical physics: "One might even say that 
the study of granular materials gives one a chance to 
reinvent statistical mechanics in a new context." 2] In- 
deed, only a few years ago the state of understanding 
of granular matter was compared to "the level of solid- 
state physics in 1930." There is no doubt that there 
are hard and deep problems associated with both types 
of systems and it may seem, at the outset, foolish to try 
to study both problems simultaneously. However, there 
have been significant advances in both fields of study 
that indicate that these problems are perhaps intimately 
related. They both deal with amorphous systems of par- 
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tides in which the dynamics is perched precariously near 
a transition between a flowing and a static state: that 
is, both systems are close to a jamming threshold where 
all dynamics ceases. One question that one can ask is 
whether there is something generic about such transitions 
so that the freezing of a liquid into a glass can profitably 
be compared to the arrest of a flowing granular material, 
or a suspension, as external stresses are reduced below 
the yield stress. In other words, can one study systems 
that can explore different states either through thermal 
fluctuations or through externally applied stresses, and 
search for unifying concepts that describe their arrested 
dynamics as different aspects of a more general "jam- 
ming" behaviorj^? 

Our approach to this problem is to describe both glassy 
systems and granular ones using the concept of a "jam- 
ming phase diagram." In such a diagram the "phase 
boundary" marks the point where the response of the 
system has become so sluggish as to make it appear solid 
on any experimental time scale. Using this framework, 
one can gain insight into the relationship between ather- 
mal jamming and thermal glass transitions and appreci- 
ate what are the control variables that govern dynamical 
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slowing down under many different conditions. In this 
paper, we describe simulations of a model liquid with 
frictionless, finite-range repulsive interactions. Because 
the potentials fall to zero at some fixed finite radius, 
such a system may be a starting point for understand- 
ing macroscopic granular or colloidal systems as well as 
liquids. For such potentials, there is one special point, 
at zero temperature and zero applied shear stress on the 
surface separating the jammed and unjammed regions, 
which has exceptional and unique properties. The goal 
of this paper is to elucidate some of the important prop- 
erties of this specific jamming transition in depth. We 
have found that the transition near this point has some 
aspects that resemble a critical point and other properties 
that are not expected for a normal second-order transi- 
tion. However, just as with a more conventional critical 
point, there is the tantalizing possibility that it may con- 
trol the region around it and thereby govern the nature 
of the entire jamming surface in the phase diagram. 

We will first describe what is meant by jamming 
and what systems may profitably be studied under this 
rubric. We will then describe the jamming phase diagram 
and show the important consequences that can be drawn 
from it. The nature of the transition at zero tempera- 
ture and zero applied shear stress will then be described 
to show why it is such an important and unique transi- 
tion. 



A. Systems that Jam 

Jamming occurs when a system develops a yield stress 
in a disordered state 0. In many cases, it is difficult 
to tell whether a system has an infinite stress relaxation 
time (and hence a yield stress), or whether it has a fi- 
nite stress relaxation time that exceeds the time scale of 
one's measurement. An alternate definition is therefore 
that jamming occurs when a system develops a stress 
relaxation time that exceeds a reasonable experimental 
time scale in a disordered state. According to these defi- 
nitions, many systems jam. Granular materials can flow 
when they are shaken or poured through a hopper, but 
jam when the shaking intensity or pouring rate is low- 
ered . Colloidal suspensions of particles are fluid but 
jam when the pressure or packing fraction is raised 
Foams and emulsions (concentrated suspensions of de- 
formablc bubbles or droplets) flow when a large shear 
stress is applied, but jam when the shear stress is lowered 
below the yield stress Q. It should be emphasized here 
that granular materials, foams and dense emulsions are 
athermal in the sense that ordinary room-temperature 
thermal fluctuations are too insignificant to allow the sys- 
tem to explore phase space. However, for other systems 
- typically those consisting of smaller particles, such as 
molecular liquids - temperature plays an important if not 
dominant role. These liquids jam (if crystallization does 
not intervene first) as temperature is lowered or density 
is increased-this is the glass transition There are a 



number of striking similarities in the phenomenology of 
these different transitions. Despite much effort, no signif- 
icant static structural signature - as opposed to a kinetic 
slowing down - of jamming has been observed experimen- 
tally in any of these systems 9J. However, we have pro- 
posed that such a signature can be observed in a quantity 
initially measured for granular materialsjldj. Another 
similarity among the different systems is that the increase 
of the stress relaxation time tends to be super- Arrhenius 
as a function of the control parameter|ll|. In addition, 
all systems show kinetic heterogeneities near the onset 
of jamming, where particle mobilities become heteroge- 
neous in space and intermittent in time |T^ . However, 
the parameters that control jamming (temperature for 
the glass transition, applied shear stress for a foam, pack- 
ing fraction for a colloidal suspension) are so different 
that it previously was difficult to see how to compare the 
jamming transitions at a quantitative level. 



B. Jamming Phase Diagram 

We proposed in Ref. 01 that different routes to kinetic 
arrest can be tied together by a "jamming phase dia- 
gram," shown schematically in Fig. ^ The shape of the 
jamming surface may be different for different systems. 
The choice of axes is dictated by the parameters that 
control the transition to jamming in the different sys- 
tems, namely temperature T, density or packing fraction 
<fi, and shear stress E. Note that T and <fi are traditional 
axes for phase diagrams, but E is not. In the unjammed 
regime, the system flows at nonzero E, so E is a non- 
equilibrium axis. Why should there be such an axis in 
the jamming phase diagram? One reason is that shear 
stress introduces fluctuations in the unjammed regime by 
forcing the system to explore different packing configura- 
tions. Recent studies show that such fluctuations can be 
described by an "effective temperature" that has many 
of the attributes of a true temperature [lj, fl5L \WL [jjj . 
Moreover, the dynamics of a sheared system whose effec- 
tive temperature is lowered toward jamming are quanti- 
tatively similar to the dynamics of an equilibrium system 
whose temperature is lowered toward the glass transition 
[IH IT^| . These results help to justify the existence of 
shear stress as an axis on the phase diagram. 

The ordinary phase diagram for the glass transition 
lies in the vertical plane coming out of the page of Fig. ^ 
namely the (1/0) — T plane. At high packing fraction 
there is a transition between a supercooled liquid and a 
glass that occurs at T g . (Although the relaxation times 
appear as if they will diverge close to the transition line, 
it is impossible in practice to track their increase past 
the times scales accessible to experiment. Thus the tran- 
sition line, T g , marks the position where the relaxation 
time has reached some large threshold. Its exact posi- 
tion may depend to a small extent on the largest time 
that an experimentalist is willing to run an experiment. 
This definition corresponds to the conventional one used 
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FIG. 1: "Jamming phase diagram" . The jammed region, near 
the origin, is enclosed by the depicted surface. The point 
labeled "J" is the boundary of the jammed region at T — 
and E = 0. Adapted from Ref. |l3j|. 



for T g in glass-forming liquids.) As the packing frac- 
tion is lowered, T g normally decreases [2(j. This glass- 
transition line is represented by the curve separating the 
jammed (i.e. glass) and unjammed (i.e. liquid) regions 
in the (l/</>) — T plane. The ordinary phase diagram 
for a foam or emulsion would be in the horizontal plane 
coming out of the page, namely the (l/</>) — S plane of 
Fig-E At fixed packing fraction, one must apply a shear 
stress higher than the yield stress in order for the sys- 
tem to flow at an experimentally measurable shear rate. 
Thus, the yield stress as a function of packing fraction 
is the curve that separates the jammed and unjammed 
regions in this plane. As the packing fraction decreases 
toward close-packing, the yield stress typically decreases, 
as indicated in Fig. ITpll 122^ . 

Mode-coupling theorists suggested years ago that the 
colloidal glass transition and molecular glass transition 
are the same despite the fact that the control variables 
are different |2^| . More recently, mode-coupling theories 
have been extended to include shear stress |24j or other 
control variables not derivable from Hamiltonians 25]. 
The jamming phase diagram suggests a reason why dif- 
ferent jamming transitions might be related, independent 
of the validity of the mode-coupling approximation. 

While it has long been recognized that temperature, 
packing fraction and stress can all control the stress relax- 
ation time, the concept of the jamming phase diagram is 
a productive way to correlate jamming in different amor- 
phous systems. The diagram implies that these three 
control parameters are important to all systems, so that 
one can study a single system as a function of all three 
variables. The diagram has proved to be a useful way to 
think about experiments, as shown recently by Trappe, 
et ai.|2^ on solidification of attractive colloids. It also 
explicitly suggests new experiments to be done. For ex- 
ample, it suggests that one should measure how the relax- 
ation time in a glass-forming liquid depends on applied 
stress. It also suggests that the introduction of a temper- 



ature to an otherwise jammed athermal system can help 
the system to flow. That is, temperature is a relevant 
variable for these transitions. This is, of course, in qual- 
itative accord with the daily experience that shaking an 
otherwise jammed material can reinitiate flow. Perhaps 
the most significant implication of the diagram is that 
the jammed region might control the behavior nearby, 
and that this is why different systems behave so simi- 
larly as they slow down on their approach to the jammed 
state. 



II. JAMMING AT POINT J 

Perhaps the most daunting problem in studying any 
jamming transition is that the jammed surface depicted 
in Fig. ^ is typically not sharp, and is defined by the sys- 
tem's relaxation time exceeding experimental time scales. 
However, there is one point on the jamming phase dia- 
gram that is well-defined [2 7| . namely, the point labeled 
"J" in Fig. n This point exists at zero temperature and 
zero applied shear stress for systems with repulsive, fric- 
tionless, finite-range potentials. This section is devoted 
to the special properties of Point J. 



A. Method 

To explore Point J, we have studied potentials of the 
following form: 

V(r ) = / 6(1 ~~ r »l ai i) a l a for ri J < ai i m 
1 13 ' \ for nj > an y ' 

where e is the characteristic energy scale of the interac- 
tion, Tij is the separation between the centers of particles 
i and j and is the sum of the radii of particles i and 
j. We study three different potentials, namely, a = 2 for 
repulsive harmonic springs, a = 3/2 for repulsive nonlin- 
ear springs that are harder than harmonic springs, and 
a = 5/2 for repulsive Hertzian interactions that are softer 
than harmonic springs [2SJ. It is important to note that 
the interactions are finite in range-particles do not in- 
teract unless they overlap. Potentials of this form were 
motivated by granular materials where particles have a 
well-defined diameter and do not interact except for a 
strong repulsive force that keeps the particles from de- 
forming too much. In our two-dimensional (2d) simu- 
lations we have used 50-50 mixtures of particles with a 
size ratio of 1.4 in order to prevent crvstalhzation [29T.l30j ]. 
The diameter of the smaller particle is denoted by a. In 
three dimensions (3d) we have studied the same bidis- 
perse mixture as well as monodisperse systems with par- 
ticle diameter a. We have studied the finite-size effects 
by varying the number of particles in the sample between 
A< N < 4096 in 2d and 3d. 

Of crucial importance is the protocol for the creation 
of configurations at T — and a given packing fraction, 
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<t>. To obtain such states, we start each simulation with a 
fixed number of particles, N, with the particle positions 
chosen completely at random (this corresponds to T = 
oo) within a square or cubic box with side length L and 
periodic boundary conditions. Starting with randomly- 
generated T — oo states guarantees that we sample all 
phase space equally. We then bring the system to the 
nearest potential-energy minimum by constantly moving 
downward on the potential energy surface. We do this 
using conjugate-gradient techniques |3l]]. Each conjugate 
gradient energy minimization is terminated when one of 
the following two stopping criteria is satisfied: 1) the 
total potential energy per particle satisfies V/N < 10~ 16 
(this corresponds to a very small pressure, p < 10 -10 ) 
or 2) V/N for successive iterations deviates by less than 
10~ 15 . This procedure brings the system extremely close 
to T = 0. Note that this procedure is identical to that 
for finding the "inherent structures" of the T — oo states 

In addition to studying the T — states generated by 
the protocol described above, we explore their properties 
by perturbing them slightly. We compress them, decom- 
press them, or apply shear strains. After each infinites- 
imal perturbation, we can again employ the conjugate- 
gradient technique. Since this technique takes the system 
to the bottom of its local potential well, the quantities we 
measure in this way are related to the static, or infinite- 
time (t = oo), response (the static bulk or shear moduli, 
Boo or Goo) of the configurations. We have also measured 
the t = moduli, Bq and Go, by measuring the response 
to a perturbation immediately after it has been applied 
(before minimizing the energy by the conjugate gradient 
technique). The shear and bulk moduli are obtained by 
measuring the response of the pressure tensor |33| 

p a fi = -L l^r, r ,-^— (2) 

to shear and compression perturbations, where Tij a is 
the a-component of Tij and d is the dimensionality of 
the system. To measure the bulk modulus we calculate 
B = <j)dp/d(j>, where the pressure is p = ^2 a p aa /d. To 
measure the shear modulus, we calculate G = dS/d'y, 
where £ = —p xy , after applying a shear shear strain in 
the ^-direction with a strain gradient in the y-direction. 
The pressure p, stress £, bulk modulus B, and shear 
modulus G are measured in units of e/<r d , lengths are 
measured in units of a, and timescales or inverse frequen- 
cies are measured in units of ay/m/e where all particles 
have equal mass m. 



B. J represents the onset of jamming for a single 
configuration 

It is important to note that each initial T = oo state 
can yield a different value of the packing fraction, cj> c , 
where the pressure and potential energy first becomes 



nonzero. Despite this ambiguity about the value of the 
threshold <f> c , we find that there are robust results when 
we measure properties as a function of <j> — 4> c , including 
scaling laws, that appear to be the same for all initial 
configurations. In Sec. Ill CI we will examine the nature 
of the distribution of these values of <p c . In this subsection 
we will show that it is possible to locate a well-defined 
onset of jamming, <p c , for each initial state. 

To test whether a given T — state is jammed or not, 
two separate criteria must be met: a jammed state must 
have a non-zero static (i.e., infinite-time) value of both 
the bulk modulus and the shear modulus. As we show 
below, for each state that we have studied, the static 
bulk and static shear moduli approach zero at the same 
density, <p c . Thus, <p c specifies the onset of jamming for 
each state. 

At T = and £ = 0, no two particles can interact if 
the density is low enough. If two particles were to over- 
lap, their repulsive potentials would simply push them 
apart during the conjugate gradient energy minimization 
process until they no longer touched. Since there is nei- 
ther thermal energy nor shear stress to compete with the 
particles' potential energy, they will never be forced back 
into contact. Thus, at sufficiently low densities there are 
no particle overlaps and the final potential energy, V, 
and the pressure, p, are both zero so that the system 
has a zero static bulk modulus. At the threshold packing 
fraction, <J) C , particles just come into unavoidable contact 
since there is no longer enough free space to allow them 
to move apart. As the system is compressed further, the 
particles overlap, the energy and pressure are nonzero 
and the bulk modulus is nonzero because the pressure 
increases upon compression. 

For each initial T = oo state, we first obtain a T = 
state using conjugate gradient minimization. For that 
T = state, we measure a precise value of <f> c , as fol- 
lows: If the configuration has zero pressure, we compress 
the system (by increasing the size of each particle by 
the same fixed fraction) in very small steps, applying 
conjugate-gradient energy minimization after each step, 
until the pressure becomes nonzero at <p c . Conversely, if 
the configuration has a non-zero pressure, we decompress 
the system in small steps, applying conjugate-gradient 
energy minimization after each step, until the pressure 
reaches zero at <f> c . We insure that the system does 
not cross over any energy barriers during these proce- 
dures by compressing (or decompressing) in successively 
smaller increments. As the density variation is made 
finer and finer, we thus make sure we end up in pre- 
cisely the same configuration for all the particles inde- 
pendent of the size of the increment. Increments were 
in the range A</> = [10 -6 , 10 -4 ], with smaller increments 
used for smaller systems and systems closer to <p c . 

At each packing fraction, we measure the static shear 
modulus, Goo, by applying a very small shear strain, min- 
imizing the energy with the conjugate gradient technique, 
and measuring the final induced stress. (Again, we insure 
that no energy barriers are crossed by applying succes- 
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FIG. 2: The infinite-time stress AE = E( 7 ) - E(0) follow- 
ing an applied shear strain 7. The resulting stress-strain 
curve is linear for sufficiently small strains and independent 
of the sign of the strain. Open (filled) symbols indicate nega- 
tive (positive) strains. These curves were generated using 3d 
monodisperse systems (TV = 512) with harmonic repulsions. 
Circles and squares represent systems with packing fractions 
cf> — <f) c — 10 -2 and 10 -4 , respectively. The solid lines have 
slopes equal to 1. The shear modulus, yield stress, and yield 
strain (where stress versus strain becomes nonlinear) tend to 
zero as 4> approaches 4> c , where (j) c is the onset of jamming for 
a given configuration. 



sively smaller increments of shear strain. The strain in- 
crements were in the range [5 x 10 -8 , 10 -5 ] with smaller 
increments used for smaller systems and systems closer 
to <j) c .) The shear modulus is calculated by measuring 
the linear relation between stress and strain, as shown in 
Fig.H 

Fig. shows the results for the pressure p as a func- 
tion of <fi — <j) c for monodisperse systems in 3-dimensions 
using both harmonic (a — 2) and Hertzian (a = 5/2) 
potentials. We also include our earlier results for bidis- 
perse systems in 2 and 3 dimensions using those same 
two potentials [23|. We find that the data for p as a func- 
tion of <fi — <fi c collapse onto a single curve for different 
initial states (each set of points corresponds to data from 
5 different states). Thus, although each initial state has 
a different value of <fi c , all states behave the same way as 
a function of <p — 4> c when compressed above </> c . 

In Fig. we show the static shear modulus, Goo, for 
the same initial states as shown for the pressure. Again, 
we find that data for different initial states collapse on 
a single curve when Goo is plotted against <p — 4> c . Note 
that <j) c was determined by where the pressure approaches 
zero, not by where the static shear modulus first ap- 
proaches zero. Thus, Figs. |21 and 01 show that the static 
shear modulus, Goo, and the pressure, p (and therefore, 
the static bulk modulus, -Boo, as well), approach zero at 
the same packing fraction, <j) c , to a precision of better 
than 2 parts in 10 for the monodisperse systems. Each 
state develops a bulk modulus and shear modulus at the 
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FIG. 3: Upper curves: Pressure p vs. (f> — <f> c for 3d monodis- 
perse (circles), 3d bidisperse (diamonds), and 2d bidisperse 
(leftward triangles) systems with harmonic repulsions (a — 
2). The solid line has slope of 2.0. Lower curves: p vs. (j> — <f> c 
for 3d monodisperse (squares), 3d bidisperse (upward trian- 
gles), and 2d bidisperse (downward triangles) systems with 
Hertzian repulsions (a = 5/2). The solid line has a slope of 
2.5. These symbols for the different systems are used through- 
out the text. TV = 1024 (TV = 512) particles were used for the 
2d (3d) systems. 



same packing fraction. This is true for all polydisper- 
sities, dimensionalities and potentials studied. Thus, (f> c 
truly marks the onset of jamming for a given initial state. 

Note that in measuring the static shear modulus, we 
apply a shear stress in a given direction. Although we 
have shown that every state studied can withstand a 
shear stress in that direction for </> > </> c , it is not obvious 
from these measurements that every state can withstand 
a shear stress in any arbitrary direction. To address 
this, we have studied the eigenvalues of the dynamical 
matrix 34] for our T = configurations with harmonic 
repulsions. We find that at least for (f> — (f> c > 1CT 6 , the 
only zero- frequency modes correspond to isolated clusters 
of "rattlers," i.e. particles that do not overlap with any 
other particles and to uniform translations of the entire 
system. The lack of any nontrivial zero-frequency modes 
shows unambiguously that the system can withstand a 
shear stress in all directions. We discuss the statistics of 
rattlers in greater detail in Section III El and the prop- 
erties of the dynamical matrix in more detail in Section 

HTTP 



C. Onset of jamming is sharp in the limit of 
infinite system size 

In the last subsection, we showed that different ini- 
tial random (T — 00) states have inherent structures 
(T = states) that jam at different threshold values, 
(j> c . Here we measure the distribution of jamming thresh- 
olds. For each system size TV and packing fraction 0, we 
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FIG. 4: Upper curves: Static shear modulus Goo vs. — 
cf>c for 3d monodisperse (circles), 3d bidisperse (diamonds) 
and 2d bidisperse (leftward triangles) systems with harmonic 
repulsions (a = 2). The solid line has a slope of 0.5. Lower 
curves: Goo vs. <j> — cj> c for 3d monodisperse (squares), 3d 
bidisperse (upward triangles), and 2d bidisperse (downward 
triangles) systems with Hertzian potentials (a — 5/2). The 
solid line has a slope of 1.0. N = 1024 (N = 512) particles 
were used for the 2d (3d) systems. 

start with at least 500 (100 for the largest system sizes) 
random (T = oo) configurations and use the conjugate 
gradient method to quench each configuration infinitely 
rapidly to T = 0. We then find the fraction of these 
final states that are "jammed", i.e., that have a finite 
pressure and static shear modulus. The resulting frac- 
tion fj of jammed states is shown as a function of <f> in 
Fig. Et a) for a two-dimensional bidisperse system and in 
Fig. EJb) for a three-dimensional monodisperse system 
with harmonic repulsions. Similar graphs were shown 
for three-dimensional bidisperse systems with harmonic 
repulsions in Ref. |27|. 

In measuring these distributions, the system remains 
at one fixed, well-defined density since we do not dilate or 
shrink the particles. Also, during the quench itself, there 
are no dynamics. The system only travels on the poten- 
tial energy surface and descends via the most rapid route 
to the nearest local potential-energy minimum. This dis- 
tribution is therefore not a function of the dynamics used 
in obtaining the final configurations but depends only on 
the fixed potential energy landscape. By starting with 
T = oo states we are sampling configuration space uni- 
formly. Thus, the result shown in Fig. EI a ) and (b) is a 
measure of the total fraction of configuration space (i.e., 
the probability) that belongs in the basins of attraction 
of final configurations that are jammed. 

Fig. EI a ) and (b) show that the fraction of jammed 
states depends sensitively on system size. For the 2d 
bidisperse system (Fig. Eta)), the curves progressively 
sharpen with increasing N, eventually approaching a ver- 
tical jump. The 3d monodisperse system (Fig. Efb)) 
shows similar behavior for N > 64. For smaller values 




FIG. 5: Fraction fj of jammed states as a function of <fi for 
(a) 2d bidisperse systems and for (b) 3d monodisperse systems 
with harmonic and Hertzian repulsions. In (a) and (b), the 
lines (downward triangles) represent potentials with a = 2 
(a — 5/2). fj for 2d bidisperse systems with a — 3/2 are 
also shown in (a) using plus symbols. Each curve represents 
a different system size N. 

of N, there is enough partial crystallization to produce 
additional structure in the curves. 

We calculate the distribution of jamming thresholds 
Pj(4> c ) by differentiating the data in Fig. El with respect 
to 4>. We find that the distributions are insensitive to 
the inter-particle potential used. This is illustrated in 
Fig. Eta) for 2d bidisperse systems at fixed system size 
N = 64. In this figure, we overlay the distributions for 
a = 5/2 (Hertzian repulsions; downward triangles) and 
a = 3/2 (plus symbols) on top of the a — 2 (harmonic; 
solid lines) distributions. In Figs. Etb)-(d), we overlay 
the distributions for a = 5/2 on top of those for a = 
2 for all systems studied (2d bidisperse, 3d bidisperse, 
and 3d monodisperse) at several system sizes N. Within 
numerical error, the different potentials yield identical 
distributions at each N. 

Fig. El also shows that it is unlikely that a jamming 
threshold (f> c will be found at very low packing frac- 
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6: (a) Distribution of jamming thresholds Pj( 
a 2d bidisperse system with N = 64 for the three differ- 
ent potentials studied (a = 3/2, 2, and 5/2). Pj{4>c) for 
(b) 2d bidisperse systems, (c) 3d bidisperse systems, and (d) 
3d monodisperse systems with harmonic and Hertzian poten- 
tials for various system sizes. In (a)-(d), the pluses, lines, 
and downward triangles represent potentials with a = 3/2, 
2, and a = 5/2, respectively. The distributions for small 3d 
monodisperse systems (N < 64) were not shown in (d) be- 
cause we wanted to emphasize the monotonic behavior of the 
peak in Pj{4>c) at large N. 



tion, where almost all states are unjammed, or at very 
high packing fraction, where almost all states are already 
jammed. For small systems, the distributions are broad; 
as N increases, they become sharper and taller. To quan- 
tify the change of the distributions with system size, we 
extract the full width at half maximum of the distribu- 
tion, w, for each N. The results are plotted in Fig. Hand 
are not monotonic in N. At very small N, there are only 
a few distinct configurations available to a static pack- 
ing, so the distribution of jamming thresholds is narrow. 
The width grows with increasing N to a maximum (near 
N = 10 for bidisperse systems and near N — 30 for 3d 
monodisperse systems). Above this value, the width de- 
creases with increasing N. At the system size where the 
distributions are widest, there is a reasonable probabil- 
ity of systems jamming at packing fractions as low as 
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FIG. 7: Width of the distribution of jamming thresholds w 
vs. the number of particles N for 2d bidisperse, 3d bidisperse, 
and 3d monodisperse systems with harmonic and Hertzian 
potentials. The solid curve has a slope of —0.55. The symbols 
have the same meaning as in Fig. ED 



roughly <j> — 0.80 in 2c? bidisperse systems and <f> = 0.58 
in 3d bidisperse and monodisperse systems. Perhaps this 
is a coincidence, but it is interesting that the value in 
3d corresponds to previous estimates of "random loose- 
packing" from experiments!^. It has been reported 
that hard particle methods (methods that strictly pro- 
hibit particle overlap) can produce jammed states with 
packing fractions that are much lower than the peak in 
the distribution of jamming onsets |36l|. However, we have 
carried out similar hard particle simulations and find that 
these low-0 states are not jammed according to our def- 
inition given above. Instead, these states are nearly un- 
jammed and fall apart when they are slightly compressed 
or sheared. 

In the large N regime, Fig.0shows that the full width 
at half maximum of the distribution scales as 



w = WqN 



-n 



(3) 



with n = 0.55 ± 0.03 and w = 0.16 ± 0.04 for all of the 
systems studied. This implies that as N diverges, the 
width approaches zero and the distribution of jamming 
thresholds approaches a 5-function. In other words, in 
the thermodynamic limit, essentially all of phase space 
jams at the same packing fraction, <j>* . This means that 
Point J in the jamming phase diagram is well-defined as 
the onset of jamming. 



D. Point J is random close-packing in an 
infinite-size system 

Our results are relevant to hard-sphere systems be- 
cause the T = configurations obtained by this proto- 
col are allowed hard-sphere configurations if none of the 
particles overlap. Thus, at sufficiently low (f>, the con- 
jugate gradient minimization technique will invariably 
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FIG. 8: Deviation of the peak in the distribution of jamming 
thresholds from its asymptotic value |0o — <j>*\ vs - L f° r 2d 
bidisperse, 3d bidisperse, and 3d monodisperse systems with 
harmonic and Hertzian potentials. The solid curve has a slope 
of \jv = 1.40. The symbols have the same meaning as in 
Fig.H 

yield allowed hard-sphere states. Our protocol yields 
special insight into the nature of random close-packing, 
a highly reproducible but heretofore somewhat vaguely- 
defined state. 

We make the connection to random close-packing by 
asking what is the limiting N — > oo value of the jamming 
threshold, <f>* . We calculate it by extrapolating the peak 
positions, 4>q, of the distributions shown in Fig. with 
respect to the system size. In Fig-IHI we plot the deviation 
of (f>o from <f)* as a function of L = N x ' d , where d is the 
dimensionality. The peak position approaches its limiting 
asymptotic value as a power law in L: 

4>o - <P* - SoL- 1 /" (4) 

By fitting to this form, we obtain v = 0.71±0.08 and So = 
0.12 ± 0.03 for all systems studied. Previously [27|, we 
obtained 0* for bidisperse systems in 2 and 3 dimensions. 
For monodisperse 3-dimensional systems, we now find 

0* = 0.639 ±0.001 (5) 

We find that <f>* does not vary with potential; this fol- 
lows from our result that the distributions of jamming 
thresholds are independent of potential (a = 3/2, a = 2, 
a = 5/2) within the uncertainty of the measurement. 
Note that the value of 4>* in Eq. 03 for monodisperse 3- 
dimensional systems is very close to what has been re- 
ported for random-close packing. Our value should be 
compared to other recent estimates of random close pack- 
ing, cj) rC p « 0.64 37, 38]. This similarity is not a coinci- 
dence. 

Random close-packing cannot be defined in a mathe- 
matically precise way because the terms "random" and 
"close- packed" are at odds with one another^. Be- 
cause the close-packing density of an fee packing is 



7r/VT8 w 0.74 > 0.64, one can always make the sys- 
tem more highly close-packed (but less random) by in- 
troducing some degree of crystalline order. How "ran- 
dom" the system should be versus how "close-packed" it 
should be is arbitrary. Torquato, et al. therefore pro- 
pose another term, "maximally random jammed state." 
By "jammed," they mean that any particle or set of par- 
ticles cannot be translated relative to any of the rest of 
the particles in the system, and by "maximally random," 
they suggest a definition based on minimization of order 
parameters characterizing the extent of crystalline order, 
bond orientational order, etc. [39]]. 

Here we suggest an alternate interpretation for ran- 
dom close-packing using the language of a "maximally 
random jammed state" but with different meanings at- 
tached to "maximally random" and "jammed." In our 
case, the value (jf is obtained by extrapolating the peak 
of the distribution of jamming thresholds to infinite sys- 
tem size. The peak of the distribution corresponds to 
the packing fraction with the maximum fraction of phase 
space (i.e., the maximum entropy) that belongs to the 
basin of attraction of jamming thresholds in the ther- 
modynamic limit. We therefore propose that another 
way to define "maximally random" is by where the en- 
tropy of initial states is a maximum, and that another 
way to define "jammed" is by the disappearance of zero- 
frequency modes of the dynamical matrix (with the ex- 
ception of isolated clusters of rattlers). This definition 
has the advantage of avoiding the order parameter de- 
scription, which will always be subject to uncertainty 
since one never knows if one has calculated the proper 
order parameter. It also provides a cleaner definition of 
the word "jammed," since it depends on nature of zero- 
frequency modes of the dynamical matrix. If one is test- 
ing whether a system is jammed by shifting particles, it 
is unlikely that one will hit on the exact combination of 
particle shifts that is characterized by the eigenvector of 
a zero-frequency mode. Finally, we note that our finding 
that virtually all initial states jam at the same value, 4>* , 
in the thermodynamic limit may explain why the value of 
random close-packing is so robust despite the fact that it 
has not been well-defined in the past. Although regions 
of the system can crystallize, such states are extremely 
rare and therefore unlikely to be observed for sufficiently 
large systems. 

The above definition of random close-packing, or the 
"maximally random jammed state," is completely well- 
defined for soft, finite-ranged repulsive potentials. What 
can be said about hard spheres? We can approach the 
hard-sphere limit by making the potential harder and 
harder-that is, by making the exponent in the potential, 
a, (see Eq. approach 0. Measuring 4>* as a function 
of a will then produce a limiting hard-sphere, value for 
random close packing. Note that our results for tp* are 
the same, within measurement error, for a = 3/2, a — 2 
(harmonic) and a = 5/2 (Hertzian). Thus, the value 
of <fi* is insensitive to a, suggesting that the hard-sphere 
limit of <f>* is the same as the value we have given in Eq.[S] 
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Of course, it is not clear that the hard-sphere limit is well- 
defined; different ways of taking the hard-sphere limit 
may lead to different results. If that is indeed the case, 
we would argue that hard spheres are a singular limit 
and thus unphysical. One should therefore concentrate 
on softer potentials for which unambiguous definitions 
can be constructed. 

Another way that has often been employed to study 
hard-sphere configurations near random close-packing 
is to conduct density ramps. For example, in the 
Lubachevsky-Stillinger algorithm |40j. a hard-sphere sys- 
tem at low packing fraction is suddenly compressed (by 
increasing the radii of all the particles at some fixed rate) 
to a higher packing fraction. In the limit of infinite 
quench rate, one finds that the system jams at a random 
close-packing density. One advantage to our protocol for 
systems with softer but still finite-ranged repulsive poten- 
tials is that, since the density is always held constant, we 
can quench the system to the final state within a fixed 
energy landscape. In the Lubachevsky-Stillinger algo- 
rithm, the energy landscape changes throughout the den- 
sity ramp because the density necessarily varies through- 
out the procedure. 

One of the strengths of our procedure is that dynamics 
has no role. If we introduce dynamics by quenching the 
temperature of the system at some finite rate, we bias 
the distributions of jamming thresholds toward higher 
values of <f>. These distributions no longer represent fea- 
tures only of the potential energy surface but now also 
depend on dynamics through the quench rate. By con- 
trast, our distributions are solely a geometric property of 
the potential energy surface. 



E. Point J is an isostatic point 

An isostatic configuration is defined by having the 
number of contacts in the system, NZ/2, equal to the 
number of force balance equations [4l| , where Z is the av- 
erage number of contacts per particle. When this occurs, 
there is a unique solution for the forces between parti- 
cles in a static packing, because the number of equations 
equals the number of unknowns. For purely repulsive, 
frictionless systems of spherical particles, the number of 
force balance equations is Nd so the isostatic condition is 
Z = 2d, where d is the dimensionality of the system. We 
find [23] that there is a discontinuous jump in Z at the 
jamming threshold, (j) c , of a given state. For cj> = <p~, 
there are no overlapping neighbors, Z = 0, while for 
4> = <t>c there are Z c overlapping neighbors. The value of 
Z c can be obtained by measuring Z at values just above 
4> c , a s shown in Fig. El The straight lines in the plots are 
fits to the data of the form 

Z-Z c = Z (<j>-<t> c f, (6) 

where £ = 0.50 ± 0.03 for all potentials, dimensions and 
polydispersities studied. 
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FIG. 9: Upper curves: Excess number of contacts per particle 
Z — Z c vs. 4> — 4> c for 3d systems: monodisperse, harmonic 
(circles); monodisperse, Hertzian (squares); bidisperse, har- 
monic (diamonds); bidisperse, Hertzian (upward triangles). 
Lower curves: Z — Z c versus </> — <j>c for 2d systems: bidis- 
perse, harmonic (leftward triangles) and bidisperse, Hertzian 
(downward triangles). N = 1024 (N = 512) particles were 
used for the 2d (3d) systems. The symbols have the same 
meaning as in Fig. 



As mentioned in Sec. Ill Bl approximately 5% of the 
particles are "rattlers" with no contacts at all, which do 
not contribute to the connected network. If we exclude 
the rattlers (so that we are only studying properties of 
the connected network) and assume £ = 0.5, then we 
obtain precise values for Z c , listed in Table [I] These 
results are consistent with Z c = 2c? in all cases, implying 
that the jamming threshold is an isostatic point. In the 
thermodynamic limit, <j) c — > (f>* , so Point J is an isostatic 
point. Note that our results for Z show that Point J is 
the only point at which the packing is isostatic; above <p* , 
we find Z > 2d so additional equations (the constitutive 
relations for the particles, which depend on the potential 
used) are needed to solve for the forces between particles. 

A more stringent condition for isostaticity is that the 
connected network (i.e. all particles in the system ex- 
cluding rattlers) has no zero- frequency modes. As dis- 
cussed in Sec. Ill Bl we have looked for zero-frequency 
modes in packings above <p c , and have tested configu- 
rations with packing fractions as little as 10~ 6 above 
<j) c . For all configurations tested, we have seen no zero- 
frequency modes except those associated with rattlers or 
with uniform translations. This suggests that Point J has 
no nontrivial zero-frequency modes. 

We have studied the fraction f r of particles that are 
rattlers as a function of <fi — <fr c for both 2d bidisperse 
and 3d monodisperse systems with harmonic interac- 
tions. We show in Figs.E3( a ) an d (b) that the fraction of 
rattlers decreases with increasing packing fraction. We 
show in Fig. HUf a) that the fraction of rattlers is inde- 
pendent of system size for N > 64 in 3d. For the 2d 
bidisperse system, we have also studied the distribution 
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FIG. 10: (a) Fraction f r of particles that are rattlers as a 
function of 4>—(f>c for a 3d monodisperse system with harmonic 
repulsions at several system sizes N. (b) f r versus (j> — 4> c for 
a, N — 1024 bidisperse system with harmonic repulsions in 
2d. (c) Number of clusters N c containing N r rattlers for five 
N — 10000 2d bidisperse systems with harmonic repulsions 
at <j> - 6 C « 10~ 2 ' 5 . 



of rattler cluster sizes. We find that most clusters have 
a single rattler and larger clusters are more rare. This is 
shown in Fig. HHT c). 



F. g(r) diverges at Point J: a vanishing length scale 

A signature of jamming at point J manifests itself in 
the pair correlation function g(r). At this point the par- 
ticles just begin to touch so an important length scale - 
the distance between nearest neighbor particles - goes to 
zero. This vanishing length scale gives rise to a diver- 
gence in g(r) in the form of (5-functions at r = , the 
sum of the radii of neighboring particles. For simplic- 
ity, we will focus on monodisperse systems. Recall from 
Sec. Ill E| that at <j>~ , there are no contacts whereas at 
Z jumps to the isostatic value Z c = 2d. This discontinu- 



ity in Z implies that there must be a (5-function in g(r) 
just at (f> c and that the area underneath this (5-function 
must be exactly the coordination number at jamming: 
Z c = 2d. This divergence is distinct from the divergence 
associated with the power-law increase above the first 
peak in g(r) (where g(r) ~ (1 — r/tj)^ 1 / 2 as r — > tr + (43) 
since that power law is integrable whereas this one has a 
nonzero area. 

Fig- Hlf al shows g(r) for a monodisperse, 3- 
dimensional system at two different values of (f>—(f> c . Note 
that as 4> approaches 4> c from above, the first peak grows 
higher and narrower. We can trace the evolution of the 
first peak by measuring its height as a function of <f> — <j) c 
(Fig- Elb)). We find that the height of the first peak at 
ro diverges as a power-law: 



g{r ) = g ((f> - (j> c )~ 



(7) 



with g = 0.90 ± 0.02 and 77 = 0.993 ± 0.002. Previous 
hard-sphere simulations 43] have measured, with much 
less precision, the height of the first peak as cf> c is ap- 
proached from below and found a similar exponent. 

In Fig. lllf c^). we plot the left-hand-width at half-height 
of the first peak of g(r) as a function oi<j)—<j) c . This width 
approaches zero as <p — > (fij as a power-law: 



A 



S = S ((j) - C ) 

where s = 0.39 ± 0.04 and A = 1.01 ± 0.005. 



(8) 



G. There is an excess low- frequency contribution 
to the density of vibrational states at Point J 

The normal modes of vibration provide a complete ba- 
sis set with which to describe the motions of the particles 
in a jammed system. There have been many studies of 
normal modes in disordered systems 0, EE EE EE EE 
EE EE EH ■ I n this section we describe the normal mode 
spectrum as a function of packing fraction above <j) c . A 
zero-frequency mode would indicate that some, possibly 
complicated, set of cooperative displacements of the par- 
ticles could be made with no cost in energy. There should 
always be d such modes corresponding to the simple uni- 
form translation of the system for each of the d dimen- 
sions. Every "rattler" particle will likewise contribute d 
zero-frequency modes. If a configuration at 4> — <j) c is 
isostatic, as we claimed in Section III El then above <f) c 
the only zero-frequency modes should be the trivial uni- 
form translations of the entire system and of the rattlers. 
As we mentioned above we have found no other, nontriv- 
ial, zero-frequency modes. On the other hand, we must 
expect some change in the nature of the low-frequency 
modes as the packing fraction for a jammed configura- 
tion is lowered toward 4> c . At that point, some extended 
mode or modes must approach zero frequency since it is 
precisely at 4> c that the system "falls apart" and becomes 
unjammed with dN zero frequency modes. How does the 
density of states evolve as <fi — 4> c approaches zero? In 
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can reduce Z all the way down to the isostatic value 
by approaching <fr+. In Fig. ^1 we show the density of 
states obtained for a monodisperse harmonic system in 3- 
dimensions with N = 1024 particles at T — 0. Fig. ITWa) 
contains the familiar result for compressed systems with 
packing fractions cf> that are far above <p c . The largest 
<t> — 4> c studied is comparable to typical liquid densities 
in an equivalent Lennard- Jones system [52j, For the very 
highest packing fractions, we see that there is an identifi- 
able region where D(u>) increases as ui 2 , as expected. As 
the packing fraction is lowered, however, we see that the 
region of lo 2 behavior shrinks, reminiscent of the results 
found in Ref. [46J. In Fig. I12f bh we show the behavior 
of D(uj) as approaches <fi c more closely. For this 1024- 
particle system, we see no sign of an ui 2 region when 
(4> — 4> c ) < 0.1. This region has presumably been pushed 
to low frequencies that are inaccessible in a system of 
this size because the excitations would have wavelengths 
that exceed the linear size of the system. Even though 
there is no cu 2 behavior at 4> — 4> c — 0.1, Fig. I12f bl shows 
that D(uS) drops as u> goes to zero. However, as <fi — <fr c 
decreases still further, this drop in D(u>) disappears. By 
<\> — (f> c — 10~ 6 there is no evidence of it at all and D(u>) 
appears to approach a constant at zero frequency. This 
striking result is unanticipated. As the packing fraction 
is lowered, the density of states approaches a limiting, 
constant, nonzero value, instead of vanishing as expected 
for long- wavelength sound modes. Thus, there is a pro- 
liferation of anomalous low-frequency modes as Point J 
is approached from above. 



FIG. 11: (a) The radial distribution function g(r) for a N = 
1024 monodisperse system with harmonic repulsions in 3d at 
4> ~ 4> c — 10 _1 and 10 -2 . The height of the first peak g(ro) 
and its left-hand-width s are defined, (b) Height of the first 
peak of g(r) as a function of 4> — 4>c f° r the same system as in 
(a). The solid line has slope —1. (c) Left-hand- width s of the 
first peak of g(r) as a function of tf> — <f> c for the same system 
as in (a). The solid line has slope 1. 



order to compute the normal modes and frequencies lo of 
the system, we diagonalize the dynamical matrix of the 
system [3lJ • The eigenvalues are the squares of the fre- 
quencies and the eigenvectors are the polarization vectors 
of the particles in each mode. 

As in a crystal one expects the low-frequency exci- 
tations to be the long-wavelength sound (longitudinal 
and transverse) modes. This assumption gives a den- 
sity of normal mode frequencies, D(lu), proportional to 
u! d ^ 1 . An earlier simulation|46| found an increase in the 
low-frequency density of states as the number of near- 
est neighbors in a glass was reduced. As we will show, 
our present results support this claim. In the previous 
study |46|, nearest-neighbor bonds were severed at ran- 
dom with some probability. Here, we control the num- 
ber of overlaps by varying the packing fraction, and we 



H. Power-law scaling near Point J 

So far, we have discussed a number of quantities that 
scale as power laws with 4> — 4> c as the jamming threshold 
is approached from the high-density side. Such quantities 
include the pressure, p (Fig. [3}, the static shear modu- 
lus, Goo (Fig. HJ and the coordination number, Z — Z c 
(Fig.®. In addition, we have shown that the width, w 
(Fig. EJ, and peak position, </> (Fig. |SJ, of the distri- 
bution of jamming thresholds display power-law scaling 
with system size. Here we discuss the power-law expo- 
nents and their implications. 

Fig-Olshows that the pressure vanishes as a power law 
as (f) — > 0+ : 

p = Pofa-&)* (9) 

The values for po and ip are listed in Table [I] Our results 
for ip are consistent with 

ip = a-l, (10) 

independent of polydispersity or dimensionality. 
The static shear modulus scales as 

G^^c^-^r (ii) 
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FIG. 12: Density of states D(u) versus u for a 3D N = 1024 
system with harmonic repulsions at packing fractions (a) far 
from <j> c and (b) close to </> c . 



where G^, and 7 are listed in Table [I] The results are 
consistent with 



7 = a- 3/2, 



(12) 



independent of polydispersity or dimensionality. 

As discussed earlier in Sec. Ill El the coordination num- 
ber Z — Z c scales as a power-law with 4> — 4>c (see Eq. |SJ) 
with an exponent consistent with 



C = 1/2, 



(13) 



independent of potential, dimensionality, and polydisper- 
sity. This result is consistent with earlier estimates from 
simulations in both 2d and 3d |27l |42j, |53J, |5J] . 

The height of the first peak of g(r) scales as a power- 
law: 



(14) 



with 77 = 0.993 ± 0.002. This was result obtained for 
a 3-dimensional monodisperse system with harmonic re- 
pulsions. Similarly, the left-hand width of the first peak 
of g(r) scales as a power-law: 



where A = 1.01 ± 0.005. 

Finally, recall the form of the fits to the width and peak 
position of the jamming threshold distributions, Eqs. 
and where the width scales asw~ N~ n and (f>* — <j)o ~ 
Fig. shows that O appears to be independent 
of potential, polydispersity and dimensionality. We find 
fl = 0.55 ± 0.03, consistent with 



1/2. 



(16) 



For the peak position, Fig.0]shows that v is independent 
of potential, polydispersity and dimensionality. We find 
v = 0.71 ± 0.08, consistent with 



2/3. 



1. Interpretation of power-law exponents 



(17) 



Some of the exponents for the scalings with cj>— 4> c are 
straightforward to understand while others are, as yet, 
without explanation. 

Pressure and bulk modulus The exponent for 
pressure, ip ps a — 1, can be explained if the system 
responds perfectly afhnely to compression. If the defor- 
mation is affine, one would expect the exponent for the 
pressure to be the same as for the force law; this argu- 
ment yields ip = a — 1. Similarly, we would expect the 
bulk modulus to behave as a power-law: 



B 



(18) 



sol 



A 



(15) 



with (3 — a — 2 because the bulk modulus is related to the 
derivative of pressure with respect to packing fraction. 
We can check to see if the response of the packing to com- 
pression is truly affine by comparing the zero-time bulk 
modulus, -Bo to the infinite-time, or static, bulk modulus, 
Boo. To obtain Bo, we apply a compression (or expan- 
sion) and measure the change of pressure without allow- 
ing any of the particles to relax their positions. By con- 
struction, the compression (expansion) is perfectly affine 
throughout the sample because we increase (decrease) 
the radii of all of the particles by the same fixed fraction. 
(This is different from how one compresses a sample in a 
laboratory experiment, where the perturbation is applied 
at the boundaries of the sample.) To obtain -Boo, on the 
other hand, we first apply the affine compression (or ex- 
pansion) , then allow the particles to shift their positions 
by minimizing the energy using the conjugate gradient 
technique. If the response to compression is perfectly 
affine, then the particles will not shift during the con- 
jugate gradient process because the energy is already a 
minimum. In that case, we would expect Bqo = Bo. The 
results are shown in Fig. For all potentials, poly- 
dispersities and dimensions studied, we consistently find 
that .Boo < -B , but that they both scale with the same 
power, consistent with f3 = a — 2. These results show that 
nonaffine deformations due to disorder in the packing do 
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FIG. 13: Zero-time (-Bo) (closed symbols) and infinite-time 
(-Boo) (open symbols) bulk moduli vs. (f> — </> c for (a) 2d 
bidisperse systems, (b) 3d bidisperse systems, and (c) 3d 
monodisperse systems with harmonic and Hertzian potentials. 
The solid curves have slopes equal to and 2.0. N — 1024 
(TV = 512) particles were used for the 2d (3d) systems. 



is the second derivative of the energy, we would expect 
the exponent for Go to be 70 = a — 2. This is indeed 
what we find, as shown in Fig. 1141 The figure shows 
that Goo < Go, as expected; the system relaxes to a 
lower value of the shear stress than it has initially. Al- 
though Lacasse, et aL|2lj| have previously pointed out 
that non-affine deformations can reduce the shear mod- 
ulus in emulsions, they did not show that the effect of 
those deformations would be to produce a power-law de- 
pendence of the shear modulus upon compression. Our 
results show that such power-law scaling exists for the 
static shear modulus and that the effect of the non-affine 
deformations is to shift the value of the exponent from 
7 = a — 2 (appropriate to the t = 0, affine situation) to 
7 rs a — 3/2 (appropriate to the t — 00 case where all 
relaxation has been allowed to take place). The effect of 
non-affine deformation is much more pronounced for the 
shear modulus than it is for the bulk modulus. In the 
latter case, the power-law exponent remained unaffected 
and only the prefactor was changed. In the case of the 
shear modulus, the non-affine deformation changes the 
scaling exponent as well as the prefactor. As the critical 
density 4> c is approached from above, the non-affine de- 
formations play a larger and larger role so that Go/Goo 
diverges at <fi c . 

Coordination number Fig.[5]shows that the coor- 
dination number scales as Z — Z c ~ (<j> — (f> c )^, where 
C is independent of potential, polydispersity and dimen- 
sionality. The fact that £ is independent of potential is 
intriguing because it suggests that £ depends only on the 
geometry of the packing. The fact that £ is also indepen- 
dent of dimensionality suggests that there is a property 
of the packing that is independent of d. 

Recent results for the pair correlation function, 
g(r), of 3-dimensional harmonic packings slightly below 
the jamming threshold show that g(r) contains a power- 
law region near r — er, where a is the sphere diameter: 



reduce the coefficient of the scaling of the bulk modulus, 
but do not change the exponent. It is not obvious why 
the exponent is unchanged. 

Shear modulus Like the bulk modulus, the shear 
modulus is also given by two derivatives of the energy. 
However, we do not find that the scaling exponent for 
the static shear modulus, 7, satisfies 7 = a — 2. Rather, 
we find 7 w a ~ 1.5 (see Eq. I12J) . To gain insight into 
this discrepancy, we have examined the zero-time shear 
modulus, Go, as well as the static or infinite-time shear 
modulus, Goo. As with the bulk modulus, to measure 
Go we first apply an affine shear strain and measure the 
resulting stress without allowing any of the particles to 
shift their positions. To measure Goo, on the other hand, 
we apply the conjugate gradient technique once the affine 
shear is applied and measure the resulting stress after 
the energy has been minimized. Since the shear modulus 



g{r) cx (1 - r/cr)- 1/2 (19) 

If one assumes an affine deformation upon compression, 
consistent with the scaling results for pressure and bulk 
modulus, then one consequence of Eq.^|is that the coor- 
dination number should increase with the power £ = 1/2, 
as we have observed. Thus, the scaling in Eq.^|is consis- 
tent with our result £ = 1/2. The origin of both results, 
however, is still not understood. 

Height and width of first peak of g(r) We find 
that the height of the first peak of g(r) diverges with an 
exponent 77 » 1 (see Ea. ll4f) and that the left-hand- width 
of the first peak vanishes with an exponent A f» 1 (see 
Eq. I15|l as <j) — > <f>+ , The fact that r\ ss A is consistent 
with our expectation that the area of the first peak is 
roughly Z c . 
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FIG. 14: Zero-time (Go) (filled symbols) and infinite-time 
(Goo) (open symbols) shear moduli vs. <f> — 4> c for (a) and (b) 
2d bidisperse systems, (c) and (d) 3d bidisperse systems, and 
(e) and (f) 3d monodisperse systems. Harmonic and Hertzian 
repulsions are labeled a = 2 and a — 5/2, respectively. In 
(a), (c), and (e) the solid curves have slopes equal to and 
0.5. In (b), (d), and (f) the solid curves have slopes equal to 
0.5 and 1.0. N = 1024 (N = 512) were used for the 2d (3d) 
systems. 



2. Discussion of finite- size scaling exponents 



We have found that there are very strong system-size 
effects. As N diverges, the width of the distribution 
of jamming thresholds vanishes as N~ n , leaving a 8- 
f unction distribution at Point J. We find that SI is very 
close to 1/2 (see Ea. lPo)) . It is not obvious that this result 
can be explained by a simple central limit theorem argu- 
ment because the packing density is a subtle property of 
the packing geometry. Independent of the explanation 
for this exponent, there are still correlations extending 
across the entire system once it is jammed. 

The peak position shifts toward the random-close- 
packing density as L~ x l v . This result suggests that there 
is a long length scale appearing in the problem near the 
onset of jamming, which scales as (<f> — 4> c )~ v ■ Note that 
our result v = 0.71±0.08 is a typical value for a correla- 
tion length exponent. 
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FIG. 15: Pressure p vs. average interparticle force (F) for a 
3d monodisperse system (N = 512) with harmonic repulsions. 
The solid line has slope equal to 1. 



I. Lack of Self- Averaging at Point J 

At Point J, there is no self-averaging in the sense that 
the average properties of a very large system are not the 
same as the average over an ensemble of many smaller 
systems at the same packing fraction. This property can 
be understood by considering a system of size N and the 
behavior as N diverges. For a finite-sized system, Fig. 
shows that there is a distribution of jamming thresholds, 
<t> c - Consider a given packing fraction, <j>, that is within 
this distribution. Some of the configurations at this <fi 
will be jammed, and others will be unjammed with p = 0. 
For an unjammed configuration, p = for every subre- 
gion of the configuration, as well. (This is exact even in 
the infinite system-size limit.) However, at the same <f> 
there will exist jammed configurations for which p > 0. 
For those configurations, we have found p > for almost 
all subregions. There are only small clusters of rattlers 
that have zero local pressures. The number of such clus- 
ters decreases rapidly with the size of the cluster (see 
Fig. UHf c)). Thus, the value of the pressure averaged over 
all configurations cannot be the same as the value of the 
pressure averaged over an arbitrary given configuration. 
As a result, there is no self-averaging. As the system 
size N increases, the distribution of jamming thresholds 
narrows. As a result, the lack of self-averaging will be ob- 
served over a smaller region of <f> that eventually narrows 
to a point (Point J) in the infinite N limit. 

The lack of self-averaging is evident in the distribution 
of inter-particle normal forces between particles, P{F) 
|27|. For a given configuration, the average interparticle 
force, (F) is directly proportional to the pressure of that 
configuration as shown in Fig. ^| for a 3d monodisperse 
system with harmonic repulsions. Depending on whether 
one normalizes the forces in a given configuration to (F) , 
the average within that configuration, and then averages 
P{F/ (F)) over many configurations, or whether one nor- 
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FIG. 16: Distribution of inter-particle normal forces for a 3d 
monodisperse (N = 1024) system with harmonic repulsions, 
(a) P(F/(F)) vs. F/(F), and (b) P(F/((F))) vs. F/«F». 



malizes the forces of all configurations to the same global 
average force ((F)), and then calculates P(F/ ((F))), one 
will get a different distribution function. This is shown 
in Fig. 1161 for a 3d monodisperse system with harmonic 
repulsions. Note that the difference between P(F/(F)) 
and P(F/((F))) is largest near = 0.636, which is near 
the peak, <f> , of the distribution of jamming thresholds 
for the 3-dimcnsional system shown (N = 1024). As 
the packing fraction is increased above </>o , the curves for 
P(F/(F)) and P(F/ ((F))) look more and more similar. 
This is consistent with the argument above, that the lack 
of self-averaging is most pronounced near the peak of the 
distribution of jamming thresholds. A simple argument 
for the shape of the tail of P(F/ ((F))) was given earlier 
271. 



J. Critical behavior near Point J 

In many ways, point J resembles a critical point. We 
have shown in Figs. 01 H El EH 113 and HI that there 
is power-law scaling near Point J of quantities such as 
the pressure, shear modulus, bulk modulus, coordination 



number, and the height and width of the first peak of 
the pair correlation function[27|]. We have also shown in 
Figs. [7| and [H] that there is finite-size scaling since the 
width and peak position of the distributions of jamming 
thresholds scale with the size of the system. This is rem- 
iniscent of behavior near an ordinary critical point. Fi- 
nally, we demonstrated in Fig. ll6l that properties such as 
the force distribution do not self- average near Point J. 
As the system size increases, the packing fraction must 
be tuned closer and closer to the peak of the distribution 
of jamming thresholds in order to see the breakdown of 
self-averaging. This is also what one expects near an 
ordinary critical point, where the temperature must be 
tuned closer and closer to the critical point as the system 
size increases in order for the correlation length to exceed 
the system size. 

The lack of self-averaging near Point J and the power- 
law scaling of the width and peak position of the jam- 
ming threshold distribution with system size all suggest 
that there is a correlation length that diverges at Point 
J. What might this length scale be? We speculate that 
there is a transverse length scale that does diverge as 
point J is approached from below. If the system is held 
at a packing fraction slightly below the critical value, 
the system is unjammed and the particles can all move 
and rearrange. However, the number of particles that 
must move in order to allow a rearrangement will de- 
pend on how close one is to the transition. Thus, in an 
infinite system, if one applies a fixed, infinitesimal ve- 
locity to a single particle we would expect the particle 
to disturb the surrounding particles as it moves. This 
disturbance will extend to a distance the transverse 
length scale, in a direction perpendicular to the applied 
force. We expect that Q will diverge as one gets close to 
the transition because as the density approaches the close 
packing value, more and more particles must rearrange 
to allow for the single particle motion in the longitudinal 
direction. The idea behind this transverse length scale is 
shown in Fig. 1171 Similar ideas are currently being ex- 
plored experimentally in granular systems with friction 
|55j and in colloidal systems [56|, |57j ■ 

One might estimate the transverse length scale by com- 
puting how many particles must move laterally in order 
to insert an extra particle. This is the parking lot model 
[58| . According to this argument, the transverse length 
scale should diverge as Q ~ ((j> c — cA) _1 ^ d_1 \ for <f> < <ft c , 
where d is the dimensionality. We note, however, that 
this result does not agree with the correlation length 
exponent that we obtained from the finite-size scaling 
analysis (see Eq. 1171 , which appears to be independent 
of dimensionality. 

Although Point J resembles a critical point, it has 
properties unlike any other critical point ever studied. 
The exponents appearing in the scaling relations are in- 
dependent of dimension but do depend on the potential. 
The former observation could be reconciled with a normal 
critical point if the upper critical dimension for jamming 
were less than 2, but then we would not expect different 
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FIG. 17: Sketch of the transverse length scale. 



potentials to yield different exponents. Likewise, if each 
different potential were in a different universality class 
and yielded different exponents, then the upper critical 
dimension should be above 3. There are other properties 
of Point J that are unusual (although not unheard of) 
for a critical point. At packing fractions below Point J, 
the pressure, shear modulus and contact number are all 
zero and the energy is zero everywhere. There are no 
fluctuations in these quantities, even infinitesimally close 
to Point J as <fi ~ * 4>* from below. In addition, there 
is a discontinuous jump in the value of the coordination 
number Z from zero to Z c at <fi* . We also note that we 
have identified a length scale that goes to zero at this 
point: the spacing between particles that form the con- 
nected network in the jammed state. This is seen in the 
divergence of the first peak of g(r) (Fig. lll|) . At a critical 
point, one expects a single divergent lengthscale and not 
a lengthscale going to zero. 

Perhaps the most disturbing feature of Point J, from 
the point of view of ordinary critical phenomena, is the 
difference in the behavior at fixed pressure and fixed vol- 
ume. At fixed volume, we observe finite-size rounding of 
power-law scaling and finite-size effects such as the lack 
of self-averaging. This is because different states have 
different jamming thresholds, cf> c . At a fixed <f>, different 
states are averaged together and the clean power-law be- 
havior we observe as a function of (f> — <j) c will be rounded. 
However, Fig. [3| shows that a fixed pressure corresponds 
to a fixed value of <fi — <j) c . When we plot quantities as 
a function of (f> — (fi c , we do not see finite-size rounding 
of power-law behavior. Even for a finite-sized system, 
the behavior of the shear modulus, etc. appears to be a 
clean power-law down to the smallest measurable values 
of cf> — <j) c . Thus, we do not see measurable finite-size ef- 
fects at fixed pressure. The divergence of g(r) also occurs 
even for a finite-sized system. These results are very dif- 
ferent from what one would expect for an ordinary critical 
point. 



III. IMPLICATIONS OF POINT J FOR THE 
GLASS TRANSITION 

We have shown that Point J marks a well-defined tran- 
sition from the unjammed to the jammed state. Because 
the conjugate gradient method allows us to probe the 
infinite-time behavior of the system, we have been able to 
show that the system develops a truly static shear mod- 
ulus at Point J. Where Point J lies with respect to the 
jamming surface depicted in Fig.^depends on one's def- 
inition. Since the glass transition line is usually defined 
as the temperature where the relaxation reaches some 
large but finite threshold value, Point J in this definition 
strictly lies within the jammed phase since the relaxation 
time there is infinite. Since Point J lies just below the 
jamming surface of the phase diagram, one might expect 
it to control behavior in its vicinity if it is indeed a critical 
point. If so, it may be the long sought-after phase tran- 
sition underlying the glass transition. In this section, we 
discuss why we suspect that the physics of point J may 
hold clues for understanding the entire jamming surface 
of Fig. ^ including the glass transition itself. 

One might wonder why Point J is important to real 
glass- forming liquids, where there are not only finite- 
ranged repulsive interactions such as those we have in- 
cluded in our calculations, but also longer-ranged attrac- 
tions. The jamming phase diagram for a real liquid would 
look quite different from the one depicted in Fig. ^ hi 
addition to the jamming surface, one has to consider the 
vapor-liquid phase coexistence curve once particles can 
attract one another. In Fig. 1181 we have sketched the 
jamming phase diagram in the T — 1/0 plane when at- 
tractions are present. For simplicity, we have explicitly 
assumed that there is no possibility of crystallization. 
(If crystallization were taken into account, then the liq- 
uid that coexists with vapor could be metastable to the 
crystal.) In Fig. ^21 the glass transition temperature de- 
creases with increasing l/cj) and eventually crosses the 
liquid- vapor coexistence region at (T x , \j$> x \ as shown. 
Once the glass transition curve crosses the left-hand side 
of the coexistence curve, which represents the lowest ac- 
cessible liquid density, a variety of states can be obtained 
depending on the quench history. The dashed part of 
the glass transition curve, which ends at Point J, is not 
necessarily accessible to systems with liquid- vapor phase 
transitions. 

Even though Point J does not necessarily exist for real 
liquids, it can still influence the glass transition. In sys- 
tems with short-ranged repulsions and longer-ranged at- 
tractions, there is still a well-defined distance at which 
the repulsion vanishes; this is the position of the mini- 
mum in the pair potential. As with the theory of liquids, 
attractions are a small perturbation to the strong repul- 
sive core; they merely hold the system at a sufficiently 
high density that the repulsions can come into play|59(. 
We therefore expect the behavior we find near Point J to 
be a good approximation to the behavior of liquids down 
to the density at which the glass transition line crosses 
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FIG. 18: A sketch of the jamming phase diagram in the 
T — l/4> plane for a system with short-ranged repulsion and 
longer-ranged attraction. For simplicity, we have assumed 
that crystallization does not occur. The jammed region lies 
underneath the curve marked T g . Beyond (T x , l/(j> x ), where 
the glass transition curve crosses the liquid density at coexis- 
tence, the available states depend on quench history. 



the liquid- vapor coexistence curve. 



Significance of divergence in pair correlation 
function 



We noted above that the first peak of the pair correla- 
tion function g{r) diverges at Point J. This has two con- 
sequences that have been observed in studies of the glass 
transition. The first has to do with the static structure 
factor, S(k), measured from scattering experiments, and 
the second has to do with the emergence of a peak in the 
distribution of normal forces, P(F), as measured exper- 
imentally in granularj^ and colloidaljilj systems, and 
numerically in previous work on models of glass-forming 
liquids^. 

At Point J, the first peak of g(r) is infinitely high and 
narrow. This property elucidates one heretofore puzzling 
aspect of studies of supercooled liquids. The static struc- 
ture factor, S(k), is related by a Fourier transform to 
g(r) so the ^-function peak in g(r) produces oscillations 
in S(k). There will not be a divergence in S(k) at any 
wavevector k. This is different from what one finds at a 
critical point where there is a diverging susceptibility at 
some value of k. (In the case of a ferromagnetic transi- 
tion, this would be the magnetic susceptibility at k = 0.) 
Thus, the signature of the transition at point J is different 
from that observed in ordinary second-order phase tran- 
sitions. As one moves away from Point J into the jammed 
region, the ((-function in g(r) broadens and decreases in 
height, but the oscillations in S(k) persist. Representa- 
tive plots of S(k) at two different values of <j> — 4>c are 
shown in Fig. ^| These are qualitatively similar to ex- 
perimental results for S(k) at high k, which also show 
oscillations jf^- This clarifies why searches for structural 
signatures of the glass transition that have examined the 




FIG. 19: Static structure factor S(k) at <j> — cj> c — 10 1 and 
10~ 4 for a 3d monodisperse system with harmonic repulsions. 



shape of S(k), either at small k or in the vicinity of the 
first peak, have not found divergent behavior. 

It has long been recognized that the first peak of g(r) 
rises and sharpens as the temperature is lowered toward 
the glass transition. However, the change of behavior 
as one crosses the glass transition is only quantitative. 
A criterion suggested many years ago(6j|, that the glass 
transition occurs when the first peak reaches a threshold 
height, seems rather arbitrary. In a previous study [Toj 
we showed that there is a qualitative change in a quantity 
closely-related to g(r). This is the distribution of normal 
forces, P(F): 



P(F)dF cx r d - l g(r)dr 



(20) 



where d is the dimensionality of the system. This quan- 
tity has been measured experimentally at the boundaries 
of static granular packings (gjj and in the interior of col- 
loidal glasses In all these studies of jammed sys- 
tems, P(F) was found to contain a peak. Our previ- 
ous studies show that a peak develops in P(F) as the 
jamming surface is approached by lowering T, increas- 
ing 0, or decreasing S 0, 0] . This signature was ob- 
served for all the potentials we have studied, including 
the full Lennard- Jones interaction, the Weeks-Chandler- 
Andersen (WCA) interaction|5°|. harmonic repulsions, 
and Hertzian repulsions. Thus, the development of a 
peak in P(F) provides a signature of the onset of jam- 
ming from purely structural data. From Eq. 1201 one can 
show that P(F) develops a peak only if the first peak of 
g(r) is sufficiently high and narrow. The criterion for a 
peak in P(F) is 



ding 1-d d 2 F/dr 2 



dr 



dF/dr 



(21) 



The fact that the onset of jamming is correlated with the 
first peak of g(r) becoming high and narrow enough sug- 
gests that the entire jamming surface may be controlled 
by Point J. 



18 




FIG. 20: Network of interparticle forces for a 2d bidisperse 
system with harmonic repulsions at <j) — <p c — 10~ 4,5 and N — 
256. The intensity of the line shading is proportional to the 
magnitude of the interparticle force. 



In order for a system to jam, it must be able to support 
shear stress for a very long time. The stress is supported 
through a network of inter-particle forces, suggesting that 
an order parameter for jamming may be found in the na- 
ture of such a network. Forces on a particle must either 
be balanced by other forces or give rise to accelerations. 
At high temperatures, there is a lot of kinetic energy and 
particles are constantly accelerated by unbalanced forces. 
At lower temperatures, however, the forces on particles 
tend to balance more because accelerations are smaller, 
and at zero temperature, forces on particles balance per- 
fectly so that the system is mechanically stable at packing 
fractions above Point J. The resulting network of forces 
at T = is shown just above the onset of jamming in 
Fig. |2SI The order parameter for the glass transition 
presumably depends on at least a three-particle quan- 
tity in order to characterize the force network. However, 
P(F), which is only a two-particle quantity, clearly cou- 
ples to the force network. A peak in P(F) reflects the 
existence of the network because the forces on all par- 
ticles can only balance if they are of roughly the same 
magnitude. This intuition highlights the importance of 
Point J to the glass transition. At T — 0, as the packing 
fraction is increased through Point J, the number of over- 
laps jumps from Z = (no force network) to Z — 2d (a 
dense force network as shown in Fig. I20|l . Thus Point J 
marks the development of a force network that supports 
shear stress. 



B. Significance of anomalous low- frequency modes 
in density of states 

Perhaps the most striking evidence that the physics at 
Point J may be related to the nature of glasses and the 
glass transition is to be found in the behavior of the den- 
sity of vibrational states at low frequencies. In contrast 
to our expectation that the density of states should vary 
as D(uj) (x lu 2 at low frequencies (in 3-dimensions), we 
find that at Point J, the density of states approaches a 
non-zero constant value as ui — > (see Sec. Ill Gil . 

We suspect that these extra low-frequency modes are 
primarily transverse in nature. It is clear from the be- 
havior of the zero- and infinite-time shear moduli, Go 
and Goo, that the transverse modes must become in- 
creasingly soft due to the relaxation allowed by non-affinc 
deformations as 4> approaches (j> c . The ratio, Go/Goo, 
diverges at this point (see Sec. I11H"H) . The bulk modu- 
lus, in contrast, does not show any particular softening 
due to non-afime relaxations and Bq/Boq, is a constant 
as <f) approaches <fi c (see Sec. IllrlH) . This suggests that 
the anomalous low-frequency modes are more transverse 
than longitudinal in character. Moreover, since Go/Goo 
diverges, and the difference between Go and Goo arises 
from spatially inhomogeneous non-affine relaxations, we 
expect that there must be significant high-wavevector 
contributions mixed into the anomalous modes |o5j . 

As <fi — > </>+, we also know that the normal modes are 
becoming more anharmonic. This was shown in Fig. [2] 
where it is clear that the linear region of the stress versus 
strain curves becomes smaller as 4> c is approached. The 
effect of this anharmonicity still needs to be determined. 

Our results, that anomalous low-frequency vibrational 
modes proliferate and herald unjamming as <f> approaches 
Point J, are of clear relevance to a large body of ex- 
perimental data on excess vibrational modes in glasses. 
Two results reflect these excess states rather directly. 
The first is the boson peak, measured by light and x- 
ray scattering |r3r| and in simulations |o7j. which indicates 
an excess of vibrational states at low frequencies, above 
those predicted by Debye behavior (D(ui) oc uJ 2 in three 
dimensions). The second is the low-temperature specific 
heat of glasses: 

c v = Arj c bycT 3 + BT + C cxcesB T (22) 

In addition to the Debye term from long-wavelength 
sound modes, there is a linear term in the specific heat 
and an excess T 3 term above that predicted by the veloc- 
ities of sound. The linear term has been ascribed to the 
existence of a new type of mode: two-level tunneling sys- 
tems P, |6^ . We note that a constant density of states 
as we have found at Point J, would, by itself, produce 
a linear term without the necessity of assuming a new 
set of tunneling excitations. However, since glasses exist 
well above Point J, we would not expect such a linear 
term to persist all the way down to zero temperature. 
Nevertheless, there is still a remarkable excess density of 
states even far away from Point J which would contribute 
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to both the excess T 3 term and to the boson peak-both 
strong signatures of glassy behavior. 



IV. CONCLUSIONS AND CONJECTURES 

We have demonstrated a number of remarkable prop- 
erties of one special point on the jamming phase diagram 
that occurs at zero temperature and zero applied shear 
stress. We have shown that this Point J is the onset 
of jamming in both the bulk and shear modulus, it is 
well defined in the N — > oo limit, and provides a clean 
definition of "random close-packing." In many ways it 
behaves as a critical point while in others it has proper- 
ties not normally associated with a second-order phase 
transition. For example, many quantities, such as pres- 
sure and shear modulus, scale as power laws with — <fr c 
but the scaling exponents depend on the inter-particle 
potential and not on dimension. 

From the finite-size scaling results that we have re- 
ported, one might also conjecture how quenched disor- 
der imposed externally (such as from pinning sites in a 
flux lattice or from optical traps in a colloidal suspen- 
sion) would affect the nature of the jamming phase di- 
agram. If we assume that the spacing between defects 
limits the correlation length in the system, instead of the 
finite size of the box that we employed in these stud- 
ies, then we would expect that the jamming threshold at 
Point J would be smeared out in much the same way as 
we find in finite-sized systems. Thus, if we were to add 
a "quenched disorder" axis on the jamming phase dia- 
gram, one of the implications of our work would be that 
as more quenched disorder is added, the distribution in 
jamming thresholds will broaden. 

Our studies here have been confined to purely friction- 
less particles. We suspect that for systems with fric- 
tional interactions, the distribution of jamming thresh- 
olds should broaden as well. This would be in accord with 
experimental observations that static frictional packings 
can exist over a wide range of densities. 

Perhaps most significant is that at Point J many of the 
properties of disordered glassy systems have their most 
pronounced expression. Just as a crystal is the most 
ordered of states, Point J may be considered to be the 
most disordered of states. As at a critical point, where 
the correlations across the entire system are most eas- 
ily observable, at Point J the nature of the disordered 
phase is most plainly seen. The constant density of low- 
frequency normal modes and the divergence in the first 
peak in g(r) are two extraordinary examples. Both are 
sensitive to global properties of the system; D(io) be- 
cause it deals with the longest wavelength modes in the 
system, and g(r) because the overlap between all par- 
ticles simultaneously goes to zero. In addition, both of 
these observations have implications for how real glassy 
systems behave. It is tempting to think that Point J may 
provide a key to understanding the nature of the entire 
surface in the jamming phase diagram and to argue that 



the properties of other glassy states should be under- 
stood as a perturbation around this "most disordered" 
of states. Thus, one might say that Point J represents 
the epitome of disorder and the essence of glassiness. 
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TABLE I: Coefficients and exponents for the power-law scaling of pressure p, shear modulus Goo, and coordination number 
Z — Z c for all systems studied. 



Power- law scaling 





System 




Quantity 








P 




Z-Z c 


d 


polydispersity 


a 


Po (±0.05) 


tp (±0.03) 


G'L (±0.05) 


7 (±0.05) 


Z (±0.5) 


C (±0.04) 


Z c (±0.02) 


2 


Bi 


2 


0.34 


1.01 


0.24 


0.47 


3.6 


0.49 


3.98 


2 


Bi 


5/2 


0.27 


1.50 


0.21 


0.99 


3.3 


0.48 


3.98 


3 


Bi 


2 


0.28 


1.03 


0.21 


0.48 


8.4 


0.47 


5.98 


3 


Bi 


5/2 


0.18 


1.51 


0.17 


1.02 


7.4 


0.49 


5.98 


3 


Mono 


2 


0.48 


1.01 


0.34 


0.49 


7.7 


0.51 


5.98 


3 


Mono 


5/2 


0.35 


1.50 


0.14 


0.95 


7.7 


0.47 


5.98 
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